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1 

Introduction 


Of the many geophysical means that can be used to probe a planet’s interior, 
seismology remains the most direct. Given that the seismic data gathered on 
the Moon over 40 years ago revolutionized our understanding of the Moon 
and are still being used today to produce new insight into the state of the 
lunar interior, it is no wonder that many future missions, both real and 
conceptual, plan to take seismometers to other planets. 

To best facilitate the return of high-quality data from these instruments, 
as well as to further our understanding of the dynamic processes that modify 
a planet’s interior, various modeling approaches are used to quantify param- 
eters such as the amount and distribution of seismicity, tidal deformation, 
and seismic structure on and of the terrestrial planets. In addition, recent 
advances in wavefield modeling have permitted a renewed look at seismic 
energy transmission and the effects of attenuation and scattering, as well as 
the presence and effect of a core, on recorded seismograms. In this chapter, 
we will review these approaches. 
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Site selection for future planetary seismology 

missions 


Typically, a minimum of three seismic stations is required for event loca- 
tion, and a fourth is required for depth determination. On Earth we take 
for granted that any given event will be relatively well-located, due to the 
comparative ease of installation of seismometers. For planetary applications, 
we cannot count on a large distribution of stations. Various factors including 
cost, difficulty of installation, instrumentation longevity, and data transmis- 
sion severely limit the number of instruments that have been or will be 
deployed on other planetary bodies. In this section, we will review various 
methods that can be employed to help determine the best landing sites for 
future planetary seismology missions, in order to maximize their scientific 
return. We focus here on the Moon and Mars, although many of these meth- 
ods are adaptable to other planetary bodies. 

Key issues that feed into the site selection of future seismometer deploy- 
ments include: 1) Theoretically establishing the seismicity and noise en- 
vironment of a planetary body, including tectonic, volcanic, artificial, and 
natural impact sources; 2) Positioning of instruments relative to regions of 
expected seismic sources to maximize the detectability of important seis- 
mic waves that traverse the interior, including body, surface, and normal 
modes; and 3) Selection of landing sites and surface materials that maxi- 
mize the longevity of a deployed instrument and minimize the seismic noise 
environment. 


2.1 The Moon 

The installation of seismometers on the Moon’s surface during the Apollo 
era provided a wealth of information that transformed our understanding 
of lunar formation and evolution. Seismic events detected by the nearside 
network were used to constrain the structure of the Moon’s crust and man- 
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tie down to a depth of about 1000 km. However, the lack of seismic ray 
paths penetrating the deepest Moon prohibited definitive identification of 
the Moon’s core. 

The presence of an attenuating region in the deepest interior, gener- 
ally interpreted as a core, has been inferred from the paucity of farside 
events (Nakamura, 2005), as well as other indirect geophysical measure- 
ments (Wieczorek et ah, 2006). In addition, more recent works have made 
progress in the recognition of core-reflected phases in the stacked Apollo 
data (Weber et ah, 2011; Garcia et ah, 2011). Such phases typically arrive 
in the coda of the main P and S arrivals, hampering their identification on 
individual seismograms. 

Because the presence of a lunar core has important implications for the 
origin and evolution of the Moon, future missions should focus on refining our 
knowledge of the deepest interior. Weber et ah (2010) devised a technique to 
map ideal landing sites for the detection of core-interacting seismic phases 
from the known distribution of lunar seismic sources. The work focused on 
reflected (PcP, ScS) and converted (PKP) phases, but is easily adaptable to 
any seismic phase. 

The method takes into account the predicted ray density and arrival am- 
plitudes (Figure 2.1), as well as the level of seismicity from the known dis- 
tribution of deep moonquakes. An example for the PKP arrival is illustrated 
in Plate 1A, using the velocity model of Nakamura (1983) modified to pos- 
sess a core with radius 340 km and P-wave velocity of 5.0 km/s. At large 
epicentral distances, PKP is predicted as a first arrival, and hence should 
be easily identifiable on future seismograms. 


2.2 Mars 

The Viking landers of the late 1970’s both carried a short-period three- 
component seismometer. On the Viking I lander, the un-caging mechanism 
failed to unlock the instrument. Viking II landed in Utopia Planitia and 
returned about 19 months of data (Lazarewicz et ah, 1981); unfortunately 
because the seismometer was not coupled to the ground, it recorded only 
wind events (Anderson et al., 1977). Since Viking, several studies aimed to 
estimate the amount of seismicity that can be expected on Mars, in order 
to optimize future seismic experiments. 

In general, three questions need to be answered to model a planets seismic 
activity: 

1. How much deformation occurs in the lithosphere? 
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Site selection for future planetary seismology missions 
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Figure 2.1 (left) Cross-section of the Moon showing the PKP ray fan for a 
source depth of 900 km. (right) Predicted arrival amplitudes for a number 
of seismic phases, including core-interacting phases. The dots show actual 
Apollo moonquake amplitudes. Note that core-reflected phases fall below 
the Apollo detection threshold. At large epicentral distances PKP is the- 
oretically detectable by the Apollo seismometers, but such source-receiver 
geometries were lacking given the limited near-site extent of the array. 


2. How does it translate into seismic events? 

3. Where do these occur? 

Philips (1991) analyzed several possible sources of deformation of the mar- 
tian lithosphere. As no plate tectonics occur today on Mars, secular cooling 
of the planetary interior is the most important source of deformation, with 
an estimated annual seismic moment release of 4.8 x 10 18 Nm, whereas tidal 
stresses can be neglected. Analyses of the accumulated slip visible on the 
martian surface today (Golombek et ah, 1992; Golombek, 2002) yield mo- 
ment release estimates of 10 17 Nm and 10 19 Nm per year. 

Knapmeyer et al. (2006) investigated not only the total amount of moment 
release, but also its distribution among individual events and the distribution 
of these events over the planet. Their model is thus not only able to link 
the seismicity to geodynamic parameters but also able to produce synthetic 
event catalogs for the evaluation of seismic experiment concepts. 

The model consists of three main parts. The annual moment release due 
to cooling of the lithosphere is first estimated. Then the moment-frequency 
relation (Gutenberg- Richter distribution) is modeled as a truncated pareto 
distribution, giving the number N of events that exceed a seismic moment 
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Mq. Finally, the seismic moment of an individual event is linked to the 
map length of the fault on which it occurred, by following the assumption 
that the stress drop is essentially constant for all quakes (e.g. Kanamori 
and Anderson, 1975). Knapmeyer et al. (2006) then mapped visible faults 
on the martian surface to distribute the individual event sizes according to 
the chosen distribution on faults that are long enough to produce them, i.e. 
small events may occur on almost any fault while large ones are restricted 
to the longest faults (Plate IB). 

Five scenarios were named to refer to the annual moment release (strong or 
weak) and the number of events (many or few) resulting from the Gutenberg- 
Richter distribution, with an average model that evenly distributed both 
factors. Translated into end-member magnitudes, the model suggests that 
the shortest fault (L = 4 km) can produce events with magnitudes up to 
Mo = 4.7, while on the longest one ( L = 1445 km), a Mq = 9.6 event is 
possible. In other words, the whole annual budget of the weak/few scenario 
may be released in a single event on each of the mapped faults, while there 
is no severe restriction for the largest possible event from the size of tectonic 
features. The recurrence times for quakes of different size, according to the 
five model scenarios, are given in Table 2.1. 

This seismicity model and the epicenter distributions derived from the 
mapped faults rely on a number of assumptions that may or may not be 
fulfilled. The amount of lithospheric contraction considered is based on a 
numerical simulation of the cooling of the lithosphere alone. Cooling of the 
entire mantle, the core, as well as a possible freezing of an inner core may 
modify the radius changes. The model for the distribution of events on the 
martian surface assumes that only the visible faults exist, and especially the 
possibility of deep sources is neglected. The latter is based on the assumption 
that deep quakes such as those that occur on subduction zones on Earth can- 
not occur in the absence of plate tectonics, and that the state of the martian 
mantle makes very deep sources similar to the lunar deep quakes impossible. 
The model by construction assumes that a strongest possible quake exists. 
The magnitude of this can of course be set to very high values. However, one 
must not confuse the annual seismic moment budget with the strongest pos- 
sible event: the largest events observed on Earth release a seismic moment 
several times higher than the average annual release. Finally, the possibility 
of seismicity linked to volcanic activity is neglected: volcanism is a mecha- 
nism to cool the planetary interior, so active volcanism would change the 
cooling rates, and related events would occur near volcanic centers. These 
two effects were ignored because the current rate of martian volcanism ap- 
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pears to be too low to be of significance for a seismic experiment in the near 
future. 

In addition to tectonic seismic events, impacts can also be used to study 
the interior of planets, and are expected to be relatively common on airless 
or thin atmospheric bodies throughout the Solar System. For the Moon, 
over 1700 impacts were recorded by the Apollo seismic array, translating 
to approximately one impact/day. Recent high resolution imaging of fresh 
craters on the surface estimates a similar cratering rate will be present on 
Mars (Daubar et ah, 2013). 
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3 


Sparse network or single-station seismology 

approaches 


Future planetary seismology missions are unlikely to have large numbers of 
instruments, as discussed previously. Therefore, many standard processing 
techniques used in terrestrial seismology will be ineffective for planetary 
applications. Specialized techniques are needed to analyze data from any 
future extraterrestrial seismic networks or single-station installations. 


3.1 Event location by P-S differential time and back-azimuth 

Future seismic missions to the Moon can take advantage of the known loca- 
tions and depths of the deep moonquake clusters detected by Apollo, which 
are presumably still active today. Future instrumentation will record seis- 
mograms from these sources on N stations, and the P and S arrival times 
t(P) and t(S) from each station will result in N x 2 unique observations. For 
a given event, we expect that with modern instrumentation we can measure 
t(S) — t(P ) with an accuracy of 5 seconds, back azimuth with an accuracy 
of 10 degrees, and the difference between the P arrivals at different stations 
with an accuracy of 1 second. With this information, we can map the clusters 
that can be uniquely identified using N desired landing sites. An example 
is shown in Figure 3.1 for N = 2. 

Most of the clusters can be uniquely identified with the S — P time dif- 
ference and back azimuth. However, the remaining clusters can possibly be 
identified by including consideration of the temporal distribution of their 
occurrence times relative to the tidal cycle. All clusters could be identified 
uniquely with four or five stations in reasonable network geometry, simply 
by location. It is important to note that the number of unique identifications 
strongly depends on the positions of the stations relative to the clusters. Fu- 
ture landing sites should be chosen carefully to take advantage of established 
seismicity. 
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Figure 3.1 Number of neighbors with which a given deep moonquake clus- 
ter can potentially be confused. A large number of clusters can be uniquely 
identified using information from only two recording stations (diamonds). 
This number increases as the number of seismic stations increases. 


3.2 Event location via timing of surface wave orbits 

With a sufficiently large event where multiple orbit surface waves are recorded, 
epicentral distance and origin time can be determined at a single station 
without any requirement for an a priori velocity model. In the traditional 
surface wave numbering scheme, orbits are numbered such that the Rayleigh 
wave propagating along the minor arc is numbered R 1, the wave propagat- 
ing along the major arc is R2, and the wave that propagates along the minor 
arc with a complete great-circle orbit is numbered R3. If we filter data con- 
taining multiple orbit surface waves using a series of narrow band filters, 
and then measure the times of the Rl, R2, and R3 peaks on the wave- 
form envelopes, we can determine the great-circle averaged group velocity 
dispersion, epicentral distance, and event origin time using simple geometry. 

For terrestrial data, we need events of approximately magnitude 6 in or- 
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Sparse network or single- station seismology approaches 


der to reliably observe the weakest arrival, R3, but this approach can be 
extended to smaller events for smaller planetary bodies. Mars, for example, 
has a radius approximately half that of the Earth, and therefore R3 would 
be expected to be about an order of magnitude larger in amplitude sim- 
ply due to geometric effects. A planet without an ocean would also likely 
have lower background seismic noise than the Earth (where the microseismic 
background noise is generated via ocean processes), and so multiple orbit 
surface waves will likely be measurable for magnitude 5 events on Mars with 
a high quality broadband seismometer, and possibly to lower magnitudes. 
This epicentral distance observation can be made independently in multiple 
frequency bands, which allows for some estimates of the error. 

Provided the surface waves are adequately recorded on the horizontal 
components, we can also determine back azimuth using the elliptical par- 
ticle motion of the Rayleigh wave by correlating rotated horizontal compo- 
nents with the Hilbert transform of the vertical component. Using such an 
approach with terrestrial data between magnitude 6 and 6.5, epicentral dis- 
tance was constrained within a degree, and origin time within 30 seconds for 
most events, despite some biases due to un-modeled 3D structure (Figure 
3.2). While such location errors would be large for the purposes of state-of- 
the art terrestrial seismology, initial tests show that they are adequate for 
performing velocity inversions using small numbers of P and S body wave 
picks, and resolve mantle velocity structure within 5%. Finally, the group 
velocity dispersion measurements included in this approach can then also 
be used to constrain shallow mantle structure using traditional surface wave 
dispersion modeling. This shallow structure then can allow better location 
of smaller events using S and P differential times as discussed previously. 


3.3 Deep structure via normal mode analysis 

For situations in which there is only a single deployed seismometer, another 
effective technique for studying deep structure is to use normal mode fre- 
quencies, which do not require knowledge of the source location. The periods 
of free oscillations of a planet are relevant to its structure. In particular the 
state of the core affects the fundamental modes of free oscillation, with the 
gravest spheroidal mode ( 0 S 2 ) being the most sensitive. 

Spectral peaks that are sensitive to mantle structure (520 mHz) can be 
detected at a noise threshold level of 10 -9 m/s 2 /Hz 1 / 2 , well within the ca- 
pabilities of modern broadband seismometers being developed for planetary 
missions (Lognonne et al., 1996). A single large quake can excite a planet’s 
normal modes, or multiple quakes with an equivalent cumulative seismic 
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Figure 3.2 Histograms of epicentral distance errors (left) and origin time 
errors (right) for 29 events between magnitude 6.0 and 6.5 recorded at 
BFO in Germany, located by picking Rl — R3 arrival times. Positive errors 
represent distances that are too large, or origin times that are too late, 
respectively (Panning et ah; 2012). 


moment can be stacked. Single-seismogram analysis requires a quake with 
seismic moment > 10 18 Nm, equivalent to a body wave magnitude mb ~ 6. 

Khan and Mosegaard (2001) attempted to investigate lunar normal modes 
using Apollo seismic data recorded from large surface impacts, in order to 
constrain seismic velocities and densities in the deep lunar interior. Their 
observations were limited to the period during which the instruments were 
operated in “flat” mode, i.e. response extended to low frequencies; only five 
events during this time were large enough. Sufficiently large events on Mars 
are similarly predicted to be rare, with less than one event per year likely 
(Golombek et ah, 1992; Knapmeyer et al., 2006); future seismic missions 
must rely on longevity to detect them. 


3.4 Location capability for future lunar seismic networks 

The traditional method for quake location converges towards the best-fit 
quake origin by minimizing the squared misfit between observed and the- 
oretical travel times. This method assumes that arrival time readings are 
Gaussian-distributed. This assumption is problematic for data from sparse 
networks; a lack of redundancy means arrival time readings with large errors 
are more likely to be included in the solution. 
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Sparse network or single- station seismology approaches 


The LOCSMITH event location scheme (Knapmeyer et al., 2008) was 
developed for sparse seismic networks on planets other than Earth. Rather 
than solving for a best-fit location, it divides the solution set into falsified and 
non-falsihed candidate locations using an adaptive grid search, and accounts 
for arrival time uncertainty using windows around the true arrival time. The 
iterative algorithm starts with an icosahedral grid over all depths of the 
planet, refining the grid size on each iteration. 

Hempel et al. (2012) applied the LOCSMITH algorithm to the Apollo seis- 
mic data. For a candidate set of grid points and detecting stations, theoreti- 
cal travel times were calculated via ray-tracing using several different known 
lunar seismic velocity models. The travel times were then compared with ob- 
served arrival time windows for both compressional wave (P) and shear wave 
phases (S). The only assumption was that the chosen time window certainly 
contains the searched arrival, which is adequate for moonquakes given that 
the first arrivals are often emergent rather than impulsive. 

If the gridded theoretical arrival times were found to lie within the arrival 
time windows for all stations and all phases, the grid point density around 
the corresponding grid point was refined for the next iteration. If not, the 
grid point was considered falsified. The solution sets of this algorithm are 
therefore clouds composed of all non-falsified grid points. The cloud size and 
shape depends on the number of observed arrivals. The best located clusters 
have ball-shaped clouds that are similar to traditional quake location error 
ellipses. Less well-located clusters have elongated clouds that resemble cones 
or bananas. The worst located clusters can only be constrained to disks that 
slice through the entire Moon. 

This correlation between cloud shape and location uncertainty can be used 
to map the detection capability of future lunar seismic networks, assuming 
the known deep moonquake clusters are still active (and we have no rea- 
son not to expect them to be). Plates 1C and ID shows the global location 
capability of an Apollo-like network, using a uniform distribution of 10,000 
synthetic events. For a relatively small-aperture near side array consisting 
of Apollo-era instruments (Plate 1C), only events occurring in the vicin- 
ity of the network will be well- located (ball-shaped uncertainty clouds). For 
a globally-distributed network consisting of modern state-of-the-art instru- 
ments (10 times more sensitive than Apollo), there is a dramatic improve- 
ment in the location capability (Plate ID). 
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Exploring the relationship between lunar seismicity 

and tidal forcing 


As discussed in Chapter 4. 5. 1.5, deep moonquakes are known to occur pe- 
riodically, with monthly (tidal) periodicity. To analyze the relationship be- 
tween lunar seismicity and tidal forcing, several previous studies created 
models of the tidal stress field within the Moon. Solid-body tides are caused 
by the gravitational influence of the Earth, resulting in the deformation of 
the Moon such that moonquakes are caused or triggered in the deep interior. 

In studies of tidal forcing, a key constraint is the focal mechanism: the fault 
parameters describing the type of failure moonquakes represent. Knowledge 
of the failure plane allows us to resolve the tidal stress into its shear and nor- 
mal components and evaluate the state of stress at the observed moonquake 
occurrence times. The normal stress resolved onto an arbitrary fault plane 
calculated at the A6 deep moonquake source is shown in Figure 4.1. A6 
moonquakes appear to occur at or near local minima in the normal stress, 
suggesting a relationship between moonquake occurrence and unloading on 
the fault. 

In terrestrial seismology, the Coulomb stress criterion predicts failure 
when a linear combination of the shear and normal stress on a fault reaches 
some critical value. Weber et al. (2009) attempted to apply a similar cri- 
terion to deep moonquakes by searching for a fault orientation such that 
moonquakes occurred when a linear combination of shear and normal stress 
on that fault reached a constant threshold. The results of their attempt for 
the Al source are shown in Figure 4.2. Fault orientations were searched in 
10° increments of dip and strike. The linear combination C of normal and 
shear stress was computed on each fault plane such that the variance of the 
ratio between C computed at moonquake times and C computed on 1-day 
intervals was minimized. The absolute minimum of this variance ratio over 
all grid spaces was selected as the preferred fault plane. 

This method worked better for some clusters than others, with a deep 
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normal stress at A6 moonquake source 
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Figure 4.1 Normal tidal stress calculated at the A6 source location, on 
a fault with a dip of 60° and a strike of 130 deg, over the course of the 
Apollo experiment. Stress is calculated at one-day intervals (line) and at 
A6 moonquake times (dots). Vertical axis unit is bars. 



a (degrees) 


Figure 4.2 Variance ratio of the linear combination of shear and normal 
stress computed at Al moonquake times to 1-day time steps, computed on 
a grid of all possible fault orientations. The best-fit fault plane occurs at 
the minimum variance ratio, and is indicated here with a star. 


minimum around the best-fit solution. For others, the definition of the best- 
fit plane was more ambiguous, with a relatively shallow minimum, or several 
nearby minima in the strike-dip grid space. Since deep moonquake failure 
planes have no surface expression, there is no way to visually verify whether 
the best-fit plane is truly representative of the physical failure mechanism. 
In addition, the quality of the lunar seismic data combined with the limited 
source/receiver geometries of the Apollo network prohibit the determina- 
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tion of deep moonquake fault parameters using first-motion polarities, as is 
typically done in terrestrial seismology. 

An alternative method for focal mechanism recovery examined the ex- 
tent to which shear (S) and compression (P) wave amplitude ratios can 
constrain moonquake fault geometry by determining whether, for a given 
cluster, there exists a focal mechanism that can produce a radiation pattern 
consistent with the amplitudes measured by the Apollo instruments (Weber 
et ah, 2012). Amplitudes were read in the ray coordinate frame, directly 
from seismograms for which the P and S arrivals were clearly identifiable 
on all long-period channels of the four Apollo stations. An empirical station 
correction accounted for site effects at the four Apollo stations, as well as 
the differences between P- and S-wave attenuation in the lunar interior. 

Instead of focusing on the best fitting solution only (as in Weber et ah, 
2009), the inverse problem was formulated using a falsification criterion: all 
source orientations that do not reproduce the observed S/P amplitude ratios 
within an error margin derived from the uncertainty of the amplitude read- 
ings were rejected. All others were accepted as possible solutions. Similarly 
to the Weber et ah (2009) method, the inversion was carried out using a grid 
search with predefined step size, encompassing all possible combinations of 
fault strike and dip, but also included slip (the direction of fault motion 
between strike and dip). 

Since the definition of a cluster implies that all events share the same 
source orientation, a combination of the inversion results of all events from 
one cluster reduces ambiguities in the inversion. Using only amplitude in- 
formation, the fault plane parameter space for a given cluster was reduced 
on average by half (Plate IE). 
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Seismic wavefield modeling 


Lunar seismograms differ greatly from terrestrial seismograms. The low at- 
tenuation of seismic energy combined with the highly scattering environment 
in the shallow subsurface results in recordings that are emergent rather than 
impulsive, with long ringing coda that obscure all but the primary seismic 
phases. For this reason, it is useful to attempt to better understand the 
transmission of seismic energy within the Moon through wavefield modeling. 
In this section we will discuss different methods for synthetic seismogram 
generation, that are useful both for comparison with existing seismograms 
(the Moon) or for the prediction of seismic signals that might be recorded 
on other planets (Mars; Europa). 


5.1 Normal mode seismology 

Complete wavefields can be predicted efficiently using summation of nor- 
mal modes calculated with a priori models of planetary structure. Such an 
approach is not well-suited for high frequency data, as the computations be- 
come more expensive and unstable for high frequencies and it is not trivial 
to include the effects of scattering and 3D structure, but they are extremely 
useful when predicting long period responses. Such calculations can be used 
in advance of particular planetary missions in order to explore likely candi- 
date seismic signals for discriminating between models of planetary interiors. 

For Europa, a moon of Jupiter, Cammarano et al. (2006) produced a 
range of possible thermodynamically consistent interior models by varying 
assumptions of initial chemistry, importance of tidal heating, and thickness 
of the ice shell overlying the likely global liquid water ocean (e.g. Kivel- 
son et al., 2000). Panning et al. (2006) used these to create normal mode 
catalogs complete up to periods of 10 seconds. By analyzing synthetic seis- 
mograms from mode summation, it was determined that long period signals 
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from internal events analogous to deep lunar events would be unlikely to be 
recorded at the surface, but fracturing events in the ice shell would likely 
produce very strong signals with characteristics strongly dependent on the 
thickness of the ice shell. In the frequency band between 1 and 100 mHz, the 
fundamental mode Rayleigh wave undergoes a transition from a nearly non- 
dispersive wave at high frequencies, to a strongly dispersed flexural wave at 
lower frequencies. This flexural wave, also observed in terrestrial sea ice at 
higher frequencies (Press and Ewing, 1951), has a characteristic dispersion 
with phase velocity proportional to the square root of frequency, and group 
velocity equal to double the phase velocity. The transition between the non- 
dispersive Rayleigh wave and the flexural mode leads to a distinctive peak 
in group velocity at a frequency which depends on the thickness of the ice 
shell. This peak means that identification of a surface wave-train, even in 
the absence of event location, should be able to determine ice shell thickness 
simply by determining the frequency of the first-arriving Rayleigh energy. 


5.2 Reflectivity and full wavefield methods 

Seismic energy transmission in the Moon is governed not only by the very low 
attenuation of waves in the surface zone, but also the intense scattering and 
reverberation of waves in the low-velocity lunar megaregolith and regolith. 
These lowered seismic velocities are similar to those observed near impact 
craters on Earth, and are related to a high amount of pore space and fracture 
density lowering the rigidity and compressibility of impact brecciated mate- 
rials (Figure 5.1a). Here we discuss methods to reproduce the long duration 
and pronounced coda observed in lunar seismic signals, and implications for 
seismic structure of planets/objects with surface megaregolith layers. 

High frequency (>lHz) synthetic wave propagation techniques are essen- 
tial for properly replicating long duration codas that relate to scattering 
and reverberation of seismic energy within the megaregolith, which normal 
mode methods don’t capture. Reverberations and trapping of seismic energy 
in a low-velocity wave channel can be achieved with a variety of 1-D wave 
propagation codes, but scattering requires a much more computationally ex- 
pensive 3-D wavefield simulation, which has only recently become available. 
Below are two examples of how high frequency wave propagation methods 
can be used to investigate the effects of low- velocity surface wave guides and 
scattering in planetary surface materials. 

Propagation of waves within the low-velocity surface materials can be 
approximated using a variety of 1-D seismic wave propagation codes. Typ- 
ically these models allow specification of a layer-cake planet and are ef- 
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Figure 5.1 Seismic structure of the lunar crust and megaregolith. A) The 
crustal velocity model after Weber et ah, (2011), showing the geology of 
the subsurface. The megaregolith extends from the surface to just above 
the structurally disturbed crust. B) 1-D wave propagation for a model with 
a 1 km thick layer of low-velocity megaregolith showing the development 
of surface waves and reverberations trapped in the layer. 


ficiently run on modern CPUs at frequencies of ~1 Hz. One example is 
GEMINI (Green’s function of the Earth by Minor Integration), made avail- 
able through the SPICE (Seismic wave Propagation and Imaging in Com- 
plex media: a European network) project. It solves the equations of motion 
for seismic source(s) and convolves the result with the Green’s coefficients 
calculated for layered 1-D media. Source functions can be either empirically 
defined or from a standard set of basic source types. The code produces seis- 
mograms for anisotropic and anelastic media for any 1-D layer-cake model 
structure; responses are calculated for a user-defined set of stations at domi- 
nant frequencies up to 10 Hz. From 1-D simulations (Figure 5.1b) we can see 
that surface waves and multiple reverberations become trapped in the low 
velocities of the lunar megaregolith, serving as a waveguide for energy trans- 
mitted into this layer. This simulates the long duration of codas but does not 
fully replicate the scattering properties observed in the lunar megaregolith. 

To properly address scattering, a high frequency, 3-D wave propagation 
code, such as WPP (Wave Propagation Program) is required (Petersson 
and Sjogreen, 2010). WPP solves the equation of motion on a 3-D finite 
difference, non-uniform Cartesian mesh, and is designed to run on parallel 
multi-CPU clusters. The user may specify elastic structures in 3D, including 
topography, layering, lateral heterogeneity, attenuation, and multiple sources 
with varying source time functions. The code is primarily used for regional 
simulations of wave propagation in complex media, the elastic response of 
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the medium is calculated for dominant frequencies up to 10 Hz or greater, 
though computational resources required grow for higher frequencies. 

In Figure 5.2, we demonstrate the effect of scattering in the megaregolith 
using WPP to obtain lunar-like seismograms. By introducing a layer of 
scattering into the low-velocity megaregolith, several lunar-like wavefield 
features, including long duration codas, partitioning of energy on all three 
components of motion, and diffuse surface waves are reproduced. Addition- 
ally, body waves entering into the coda layer from below become scattered 
at the surface only. As computational resources grow, fully global simula- 
tions at high frequency are becoming available, allowing us to investigate 
the global effects of scattering at very long distances. Scattering effects will 
vary from world to world, with airless bodies with well-developed megare- 
goliths (Moon, Mercury, asteroids) having the highest proportion of scatter- 
ing within their crusts relative to bodies with atmospheres and intragranular 
fluids and ices (Mars, Earth, Venus, the outer satellites). WPP is appropriate 
for any planetary object with a solid surface and can be used to approximate 
the scattering from rocky and icy worlds. 
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Figure 5.2 Seismic scattering in the lunar crust and megaregolith for an 
impact source. A) Cross section through a lunar model with 25%, 200 
m scale-length, Von Karman random heterogeneity distribution (Frankel 
and Clayton, 1986) within a 3 km thick megaregolith and time snapshot 
of vertical displacement at 40 seconds into a WPP simulation. Positive 
amplitudes are black, negative white. B) Seismograms for a model with no 
scattering in the subsurface (left) and 25% heterogeneity (right) at 50 km 
from the source (location given in part A). 






Plate 1 . A) Nearside map of PKP ray density from the known distribution of deep moonquake clusters as a 
function of landing site coordinates, taking predicted cluster occurrence activity and arrival amplitudes into 
account, shown in logarithmic scale. For this particular phase, landing sites near the limb are favored, particu- 
larly in the southwestern quadrant of the Moon, where the likelihood of detecting PKP from the northeastern 
farside events is greatest. B) Map showing the probability that a M w <4.8 quake is assigned to 1 km 2 of 
surface, plotted on a 2.5° x 2.5° grid (Knapmeyer et ah, 2006). C) Comparison between synthetically derived 
location clouds and the Apollo data, for an Apollo-like network. Dots represent ball-shaped located clusters 
of deep moonquakes. Triangles are station locations (Hempel et ah, 2012). D) Location capability of a 
network with four equally distributed, state-of-the-art seismometers (Hempel et ah, 2012). E) Least-squares 
misfit between observed and calculated amplitudes on a regular grid of fault orientations averaged over all 
useable events from the A6 source region. The accepted regions (green) are selected by highlighting fault 
orientations which predict amplitudes within two standard deviations of the measured values. 
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